function [X,Y,Z] = computeMeshforPlot(robot_input,field_limits,dim)
% compute the full distribution of the objective for all the ground points.
    if size(robot_input,2) == 1
        robot = reshape(robot_input,numel(robot_input)/3,3);
    else
        robot = robot_input;
    end
    h = robot(:,end); 
    r = coverageModel(h,20,'reciprocal');
    dis_num = dim;r_num = size(robot,1);
    x_space = linspace(field_limits(1),field_limits(2),dis_num);
    y_space = linspace(field_limits(3),field_limits(4),dis_num);
    [X,Y] = meshgrid(x_space,y_space);
    points_xy = [reshape(X,dis_num^2,1) reshape(Y,dis_num^2,1)];
    count_num = zeros(dis_num^2,r_num);
    for i = 1:r_num
        posi = points_xy - robot(i,1:2);
        dis = vecnorm(posi,2,2);
        count_num(:,i) = (dis< r(i));    
    end
    
    z = zeros(dis_num^2,1);
    r1 = sum(count_num,2) == 1;  z(r1) = 0.05;
    r1 = sum(count_num,2) == 2;  z(r1) = 0.1;   
    r1 = sum(count_num,2) >=3; r3 = kron(r1,ones(r_num,1)); 
    r2 = reshape(count_num',numel(count_num),1);  r2 = r2&r3;    
    robot_rel = kron(ones(dis_num^2,1),robot);
    points = [kron(points_xy,ones(r_num,1)) zeros(r_num*dis_num^2,1)];

    G = (points - robot_rel);
    G = G.*(r2)./vecnorm(G,2,2);   
    
    k1 = find(r1==1);
    for j = 1:size(k1)
        i = k1(j);
        Gi = G(r_num*(i-1)+1:r_num*(i-1)+r_num,:);
        z(i) = sigmoid(  pdop(Gi)  );
    end
	Z = reshape(z,dis_num,dis_num);
end

